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ABSTRACT 

We introduce the Marenostrum-MultiDark Simulations of galaxy Clusters (MUSIC) dataset. 
It constitutes one of the largest sample of hydrodynamically simulated galaxy clusters with 
more than 500 clusters and 2000 groups. The objects have been selected from two large 
N-body simulations and have been resimulated at high resolution using Smoothed Particle 
Hydrodynamics (SPH) together with relevant physical processes that include cooling, UV 
photoionization, star formation and different feedback processes associated to Supernovae 
explosions. 

In this first paper we focus on the analysis of the baryon content (gas and star) of clusters 
in the MUSIC dataset both as a function of aperture radius and redshift. The results from our 
simulations are compared with a compilation of the most recent observational estimates of the 
gas fraction in galaxy clusters at different overdensity radii. We confirm, as in previous sim- 
ulations, that the gas fraction is overestimated if radiative physics is not properly taken into 
account. On the other hand, when the effects of cooling and stellar feedbacks are included, the 
MUSIC clusters show a good agreement with the most recent observed gas fractions quoted 
in the literature. A clear dependence of the gas fractions with the total cluster mass is also 
evident. However, we do not find a significant evolution with redshift of the gas fractions at 
aperture radius corresponding to overdensities smaller than 1500 with respect to critical den- 
sity. At smaller radii, the gas fraction do exhibit a decrease with redshift that is related the 
gas depletion due to star formation in the central region of the clusters. The impact of the 
aperture radius choice, when comparing integrated quantities at different redshifts, is tested. 
The standard, widely used definition of radius at a fixed overdensity with respect to critical 
density is compared with a definition of aperture radius based on the redshift dependent over- 
density with respect to background matter density: we show that the latter definition is more 
successful in probing the same fraction of the virial radius at different redshifts, providing 
a more reliable derivation of the time evolution of integrated quantities. We also present in 
this paper a detailed analysis of the scaling relations of the thermal SZ (Sunyaev Zel'dovich) 
Effect derived from MUSIC clusters. The integrated SZ brightness, Y, is related to the cluster 
total mass, M, as well as, the M — Y counterpart which is more suitable for observational 
applications. Both laws are consistent with predictions from the self-similar model, showing a 
very low scatter which is ai og y — 0.04 and even a smaller one {<J\ og m — 0.03) for the inverse 
M — Y relation. The effects of the gas fraction on the Y — M scaling relation is also studied. 
At high overdensities, the dispersion of the gas fractions introduces non negligible deviation 
from self-similarity, which is directly related to the f gas — M relation. The presence of a 
possible redshift dependence on the Y — M scaling relation is also explored. No significant 
evolution of the SZ relations is found at lower overdensities, regardless of the definition of 
overdensity used. 

Key words: methods: iV-body simulations - galaxies: clusters - galaxies: scaling relations - 
cosmology: theory - Sunyaev-Zel'dovich effect 
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1 INTRODUCTION 

Galaxy clusters are the biggest gravitationally bound objects of 
the Universe and constitute one of the best cosmological probes 
to measure the total matter content of the Universe. However, the 
total mass of these objects cannot be directly measured. It must be 
inferred from other observational quantities (X-ray or SZ Surface 
Brightness, Lensing distortions or number of galaxies). In all cases, 
one has to relate these quantities with the total mass of the sys- 
tem. Due to the complex physics involved in the processes of clus- 
ter formation, hydrodynamical numerical simulations have been a 
fundamental tool to calibrate mass proxies, define new ones, and to 
study the systematics affecting observational measurements. They 
are also indispensable to deeply study the formation and evolution 
of clusters of galaxies and all their gas-dynamical effects ( |Borgani| 
|& Kravtsov|2009} . 

The big progresses achieved in the last years by numerical 
simulations are well represented by their use to describe and to 
study X-ray temperatures and their relation with cluster gas mass 
( |Ettori et al.|2004||Muanwong et al.|2006||Nagai et al.|2007a| > as 
well as to compare numerical predictions with observed tempera- 
tures ( |Loken et al.|2002||Borgani et al.|2004||Leccardi & Molendi| 
|2008l > , gas profil es jRoiicarelli et al.|2006||Croston et al.|2008) , or 
pressure profiles ( |Arnaud et al.|2010 . 

In an ideal scenario one would need to have a large sample of 
simulated galaxy clusters with enough numerical resolution (both 
in mass and in the gravity and pressure forces) to accurately re- 
solve the internal substructures and with a detailed modelling of 
the most relevant physical processes. The best way to achieve this 
goal would be by simulating large cosmological boxes ([Borgani 



et al.||2004| jGottlober & Yepes||2007l |Burns et al. 2008; Hartley 



et al.||2008[ |Boy lan-Kolchi n et al.||2009fr . Unfortunately, due to 



the large computational demand of these simulations, one needs 
to find a compromise between the three main ingredients: volume 
size, mass resolution and physical processes included . A possi- 
ble solution to the computational problems related with scalability 
of the present-day hydro codes is to proceed, mimicking the ob- 
servations, by creating a catalogue of resimulated galaxy clusters 
that are extracted from low resolution N-body simulations. The re- 
gions containing clusters of galaxies are then resimulated with very 
high resolution, adding only gas physics in the resimulated areas 



(Puch wein et al.||2008|Dolag et al.|[2009) |Lau et aLp0 09 Fabjan 
|et al.|2010) . This so-called 'zooming' technique permits to simu- 
late thousands of clusters basically independently from each other 
with less computational cost than a full box hydrodynamical sim- 
ulation of the same resolution. By selecting all the objects formed 
in a given volume above a given mass threshold, mock volume lim- 
ited sample catalogues can be generated and used in the study of the 
properties and interrelations of the different scaling laws of galaxy 
clusters. Following this procedure, we have generated the MUSIC 
(Marenostrum-MultiDark Simulation of galaxy Clusters) dataset, 
a large sample of simulated clusters of galaxies, composed of ob- 
jects extracted from two large box cosmological simulations: the 
MareNostrum Universe and the MultiDark simulation. We selected 
all the clusters using criteria based on mass (selecting all the clus- 
ters having a total mass larger than 10 15 /i _1 M ) or on morphol- 
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ogy (selecting groups of clusters corresponding to different mor- 
phology classes, bullet-like clusters and relaxed clusters). All these 
objects have been resimulated with SPH particles, radiative phys- 
ical processes and star formation prescriptions, improving by an 
order of magnitude the resolution with respect to the original sim- 
ulations, as described in the next section. From this database, we 
will obtain mock observations for X-rays, SZ, lensing as well as 
optical galaxy counts. This will allow us to study the interrelations 
between the scaling laws associated to the different observables. In 
this paper we base our attention on the properties of the baryon con- 
tent and the SZ effects and will leave a more detailed analyzes of 
the relations with the other observed properties for a further work. 

In the last few decades the SZ effect ( Sun yaev & Zeldovich| 
|1970[ > has become one of the most powerful cosmological tool to 
study clusters of galaxies, as well as the nature of the dark matter 
and dark energy components of the Universe. The physical pro- 
cess of the SZ effect is the diffusion of CMB (Cosmic Microwave 
Background) photons with a hot plasma due to inverse Compton 
scattering. The thermal component of the SZ effect is largely en- 
hanced by the presence of clusters of galaxies, the most massive 
bound objects in the Universe, where plasma is in hydrostatic equi- 
librium inside the gravitational potentials of dark matter. The Intra 
Cluster Medium (ICM) composed by high energy electrons consti- 
tutes an ideal laboratory to investigate the SZ effect. The brightness 
of the SZ effect turns out to be independent of the diffuser position, 
thus making it the best tool to find galaxy clusters at high redshift. 
Moreover, the SZ flux collected from the cluster region is propor- 
tional to the total thermal energy content, with a weak dependence 
on the complex physical processes acting at the inner regions (e.g. 
cooling flows, galaxy feedbacks etc) which mostly affect the X- 
ray luminosity. Therefore, these two properties: redshift indepen- 
dence and low scatter mass proxy makes the integrated Compton 
y-parameter an efficient high-z mass-estimator. 

Under the hypothesis that the evolution of galaxy clusters is 
driven mainly by gravitational processes (Kaiser 1986) and assum- 
ing hydrostatic equilibrium and an isothermal distribution of dark 
matter and ICM (Bryan & Norman 1998) it is possible to derive 
simple scaling power-laws linking cluster properties: the so-called 
self-similar scaling relations. In the case of SZ science, the relation 
linking the SZ brightness with the cluster total mass, the Y — M 
scaling law, is continuously under analysis to test for its robustness, 
allowing the application of the SZ effect as a mass-finder. 

Observational studies started to collect data of a few clus- 
ters, mostly those with high X-ray luminosities (and therefore high 
masses) ( Benson et al. 2004; Morandi et al. 2007 Bonamente et al. 
|2008||Vildilinin et al.|2 009; Arn aud et al.|2010J . Recent large sur- 
veys have shown the possibility of detecting undiscovered clusters 
only through SZ effect observations, as in the claims by South Pole 
Telescope (SPT,|Staniszewski et al.|2009[|Vanderlinde et al.|2010| 
l |Andersson & SPT Collaboration|2010) ; |Williamson et al.|2011| >, 
Atacama Cosmology Telescope (ACT, Marriage et al.||201 1| |Se-| 
and Planc k (|Pla nck Collaboration et al.||2011a| 
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Planck Collaboration et al. 2011c; Planck Collaboration et al. 
|201 lb |Planck Co llaboration et al. 201 Id). The possibilities to ex- 
plore more distant objects or to deeply map single cluster morphol- 
ogy are planned with the ongoing higher angular resolution projects 
like AMI jZwart et al.|2008] >, AMiBA ( |Lo et al.|200l) , MUSTANG 
l Korngut et al.pOTT] ), OCRA ( |Browne et al.|2000^ CARMA pi>] 
|chovej et al.|2012 1 or the in-coming projects like the ground based 
C-CAT ( [Golwala et al.|2007|, the upg raded with new spectroscopic 
capabilities MITO l |De Petris et al.||2007) and the balloon-borne 
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OLIMPO jMasi et al.|2008| > or the proposed satellite mission Mil- 
limetron (http : / / www . sron . rug . nl/millimetron). 

More recent blind-surveys carried out by SPT, ACT and 
Planck enlarged the existing dataset confirming the self-similarity 
in the sample at least for the massive clusters (Marriage et al. 201 1 



Ande rssonetal,|2011[ |Planck Co llaboratio n et al.|201 lb) . Numer 
ical simulations have shown that Y (the integral of the Compton y- 
parameter over the solid angle of the cluster) is a good mass proxy 
qda Silva et al.|2004| |Motl et al.|2005]|Aghanim et al.|2009| > and 
that the slope and the evolution of SZ scaling relations are appar- 
ently not affected by redshift evolution and cluster physics. In order 
to find an X-ray equivalent of the SZ integrated flux, the numerical 
simulations also led to the introduction of a new mass proxy: the 
Yx parameter, defined as the product of the cluster gas mass and 
its tempearture (Kravtsov et al. 2006). 

Many works based on large A r -body cosmological simulations 
have already studied the impact of gas physics on SZ scaling rela- 
tions (da Silva et al.|2000| [White et al.|2002l|da Silva et al.|2004 



Motl et al. 2005; Nagai 2006, Bonaldi et al. 2007; Hallman et al 



20071 |Aghanim et al.|2 009 ; Battaglia et al.|201 l[|Kay et al.|2012> , 



such as the effects introduced by clusters with disturbed morpholo- 
gies JPoole et al.|2007||Wik et al.|2008l |Yang et al.|20T0} |Krause] 
|et al.|2012| , showing that the self-similar model is valid at least up 
to cluster scales (even if with some differences introduced by the 
different models used in the simulations). 

In what follows, we will present an analysis of the relation 
between the SZ integrated flux and the total mass, Y-M, in MU- 
SIC clusters and will compare it with the predictions of the self- 
similar model. We also study the possible biases introduced by the 
common assumption of considering quantities (M or Y) integrated 
inside a radius defined by a fixed overdensity with respect to the 
critical density instead of a more suitable definition with respect to 
the background density whose value depends on redshift. 

We will focus our analysis of the Y-M scaling relation on the 
most massive objects of the MUSIC dataset which constitutes an 
almost complete volume limited sample. Therefore, only clusters 
with virial masses M v >5x 10 fo Mq, are considered in this 
paper. We will extend our analysis to a wider range of masses in an 
upcoming work. 

The paper is organized as follows. In Section 2 the MUSIC 
database is described. The baryon content in the clusters at differ- 
ent overdensities is presented in Section 3, together with a study 
of its evolution with redshift and a comparison of our numerical 
results with the most recent observational estimates. In Section 4 
the Y-M scaling relation is computed from the MUSIC clusters 
and compared with the predictions from the self-similar model. We 
also discuss the validity of the integrated Y as a proxy for the clus- 
ter mass. In Section 5 we focus on the impact of the gas fraction 
on the Y-M relation together with the dependence of the gas frac- 
tion on the cluster mass. Finally in Section 6 the redshift evolution 
of the Y-M relation is discussed. In Section 7 we summarize and 
discuss our main results. 



2 THE MUSIC DATASET 

The MUSIC project consists of two sets of resimulated clusters ex- 
tracted from two large volume simulations: 

• The MareNostrum Universe, a non-radiative SPH simulation 
with 2 billion particles (2x 1024 3 gas and dark matter ) in a 500/i _1 
Mpc cubic box ( |Gottl6ber & Yepes|2007l >. 



• The MultiDark Simulation, a dark-matter only N -body simu- 
lation with 2048 3 particles in a 1 h~ 1 Gpc cubic box iPrada et al. 
[20T7) . 

These two simulations have slightly different cosmologies. 
The Marenostrum Universe (MU) was made with the cosmo- 
logical parameters that were compatible with WMAP1 results ( 
f2 M =0.3, ^6=0.045, n A =0.7, <78=0.9, n=1.0, h=0J) while the 
MultiDark run was done using the best-fit cosmological parame- 
ters to WMAP7 + BAO + SNI ( D. M =0.21, n 6 =0.0469, «a=0.73, 
(78=0.82, n=0.95, h=0J) ( |Komatsu et al.|2011| >. 

The procedure to select the interested objects in those two 
simulations was also different. For the MareNostrum clusters 
we mainly selected them based on the dynamical state (relaxed 
vs bullet-like clusters). For the MultiDark clusters we made a 
mass limited selection, taking all clusters with masses above 10 15 
fr -1 M Q at2 = 0. 

The zooming technique (Klypin et al. 2001 1 was used to pro- 
duce the initial conditions of the resimulated objects. Namely, we 
first found all particles within a sphere of 6 Mpc radius around the 
center of each selected object at z=0 from a low resolution ver- 
sion (256 3 particles) of the two simulations. This set of particles 
was then mapped back to the initial conditions to find out the La- 
grangian region corresponding to a 6 h~ 1 Mpc radius sphere cen- 
tred at the cluster centre of mass at z = 0. The initial conditions 
of the original simulations were generated in a finer mesh of 2048 3 
(for the MareNostrum) and 4096 3 (for the MultiDark) sizes. There- 
fore, we could improve the mass resolution of the resimulated ob- 
jects by a factor of 8 with respect to the original simulations. We 
kept the highest mass-refinement level within the Lagrangian re- 
gion of each cluster and then cover it with shells of increasing mass 
particles down to the lower resolution level of 256 3 . To avoid prob- 
lems with periodic conditions, we always recentre the simulations, 
each resimulated cluster always located at the center of the cor- 
responding box. Thus, for the MultiDark clusters, we have dark 
matter particles of 5 different mass refinements ( from 4096 3 to 
256 3 ) while for the MareNostrum clusters we have 4 different mass 
species. The gas SPH particles were added only to the highest re- 
finement level. The SPH particle positions were slightly displaced 
from their parent dark matter by 0.4 times the mean inter particle 
distance in the 3 spatial directions and they were given the same 
initial velocity as their dark matter counterparts. 

The parallel TreePM+SPH GADGET code ( |Springel|[2005] > 
was used to run all the resimulations. We accounted for the ef- 
fects of radiative cooling, UV photoionization, star formation and 
supernova feedback, including the effects of strong winds from su- 
pernova, in the same way as described in Springel & Hernquist| 
(2003 ) model. We refer the reader to this reference for a detailed 
description of the model. Here we simply provide the values corre- 
sponding to the typical parameters of the model. Stars are assumed 
to form from cold gas clouds on a characteristic timescale t* and 
a certain stellar mass fraction f3 is instantaneously released due to 
supernovae from massive stars (M > 8 Mq). We adopt a fraction 
of /?=0. 1 which is consistent with assuming an Universal Salpeter 
IMF with a slope of -1.35 in the limits of 0.1 M and 40 M Q . 
In addition to this mass injection of hot gas, thermal energy is also 
released to the ISM by the supernovae. This feedback energy heats 
the ambient gas and cold clouds can be evaporated inside hot bub- 
bles of supernovae. All these assumptions lead to a self-regulating 
star formation, where the growth of cold clouds is balanced by su- 
pernova feedback of the formed stars. The number of collision-less 
star particles spawned from one SPH parent gas particle is fixed to 
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2 , which means that each stellar particle will get half of the orig- 
inal parent SPH mass. In addition to the thermal feedback, some 
sort of kinetic feedback is also accounted for. Gas mass losses due 
to galactic winds, M w , is assumed to be proportional to the star 
formation rate Msfr, i.e. M^^Msfb, with r\=2. Additionally, the 
wind contains a fixed fraction e of the total supernova energy, which 
is assumed to be e = 0.5. SPH particles near the star formation re- 
gion will be subject to enter in the wind in an stochastic way. Those 
particles eligible for winds will be given an isotropic velocity kick 
of w=400 km/s and will freely travel without feeling pressure forces 
up to 20 kpc distance from their original positions. 

Using the aforementioned techniques we have been able to 
resimulate a total of 700 different galaxy clusters from the two large 
box simulations. This constitutes one of the largest samples of res- 
imulated cluster-sized objects. We used two different selection cri- 
teria, based on dynamical state and on mass cuts in order to explore 
scaling relations dependencies on different aspects (mass, morphol- 
ogy, redshift evolution, overdensity). All these clusters were also 
resimulated with non-radiative physics, using only SPH and grav- 
ity forces, so we can also account for the effects of cooling, star 
formation and supemovae feedbacks on the observed properties of 
the clusters. 

The database of MUSIC clusters is publicly avail- 
able in a SQL query format through the website: 
|http : / /music . ft . uam . es| The database will compre- 
hend also the initial conditions of all MUSIC objects, in order 
to give the possibility to resimulate the clusters with other hydro 
codes and/or with different modelling for the radiative processes 
(e.g. AGN feedback, MHD, Cosmic ray pressure, etc). 

In what follows we describe in more detail the two main sub- 
samples in which the MUSIC database is divided. 

2.1 MUSIC-1: MareNostrum resimulated clusters 

This first subset is composed of 164 resimulated clusters extracted 
from the MU. 

There are more than 4000 clusters of galaxies (objects larger 
than 10 14 h~ 1 M e ) in the whole MU simulation. Out of this 
large number of clusters, we selected 82 'bullet-like' clusters. The 
1ES0657-556 cluster (Clowe et al. 2006), generally known as the 
bullet cluster, is considered one of the best astrophysical tools to 
study the nature of dark matter. This cluster, located at z — 0.296, 
is composed by a smaller merger which has just crossed a bigger 
object and is now moving away. The core of the plasma cloud of the 
merging sub-cluster, 'the bullet', has survived the passage. Due to 
the different dynamics of the collision less dark matter component 
and the x-ray emitting gas, it is expected that the two components 
have a spatial offset. By combining lensing and X-ray observations, 
this behaviour is seen explicitly in the bullet cluster. The relative ve- 
locity of the two merging systems is a quantity that is not directly 
measured. It is derived from the shape of the "bow shock" clearly 
seen in X-ray images, using the results from different hydro simu- 
lations of merger clusters ( Springel & Farrar 2007 , Mastropietro & 
|Burkert|2008l >. But a more direct measure of the effects of a cluster 
merger, in the gas and dark matter components, is the projected 2D 
separation between the peak of the two distributions (gas and dark 
matter). Forero-Romero et al. (2010) have used this measurement 
to characterize the number of bullet-like clusters in the MU. They 
define a "bullet" cluster if the 2D separation between the peak of 
dark matter and the peak of gas distribution is at least 200 kpc (as 
in the real bullet cluster). 

In order to have a morphological counterpart to these ex- 



tremely disturbed systems, we also selected another 82 clusters 
which exhibit the following relaxed conditions: the displacement 
between dark and baryonic matter is smaller than 200 kpc and they 
are composed of a single massive cluster with no big substruc- 
tures (i.e. all substructures inside the virial radius must have masses 
smaller than 10 per cent of the total mass). 

The relaxed clusters were chosen to have masses similar to 
those of the bullet systems: in this way every bullet-like cluster 
has a relaxed companion of the same mass. All the selected clus- 
ters were resimulated using the zooming technique and radiative 
physics, including radiative cooling, heating processes of the gas 
arisen from a UV background, star formation and supernovae feed- 
back explained above. The particle mass of dark matter was set 
to moM=l-03x lO 9 /i _1 M and the mass of SPH gas particle to 
m 9as =1.82x 1O 8 /i _1 M . Due to the relative small size of the com- 
putational box of the MU simulation, very few clusters of MUSIC- 

1 subset have a mass at the virial radius bigger than 10 15 h~ M . 
The cluster masses of MUSIC-1 span from M v =2x lO 15 /i _1 M to 
M I , = lxlO 14 ft" 1 M . 

2.2 MUSIC-2 : MultiDark resimulated clusters 

The second, and more numerous, subset of MUSIC consists of 
a mass limited sample of resimulated clusters selected from the 
MultiDark Simulation (MD). This simulation is dark-matter only 
and contains 2048 3 (almost 9 billion) particles in a 
cube. It was performed in 2010 using ART (Kravtsov et al. 
|1997| > at the NASA Ames Research centre. All the data of this 
simulation are accessible from the online MultiDark Database 
(www.MultiDark.org). We have selected, using a low resolution 
(256 3 particles) version of the MD simulation, all the objects more 
massive than 10 15 Ii~ 1 Mq. In total, we found 282 objects above 
this mass limit. All these massive clusters were resimulated both 
with and without radiative physics using the techniques explained 
above. The mass resolution for these simulations corresponds 
to m Dil /=9.01xl0 8 /i~ 1 M Q and to m S PH=1.9x lO 8 ft _1 M . The 
gravitational softening was set to 6 ft -1 kpc for the SPH and dark 
matter particles in the high-resolution areas. Several low mass clus- 
ters have been found close to the large ones and not overlapping 
with them. Thus, the total number of resimulated objects is con- 
siderably larger. In total we obtained 535 clusters with M>10 14 
ft _1 M at z = and more than 2000 group-like objects with masses 
in the range 10 13 /i _1 M Q <M„ <10 14 /i _1 M . All of these res- 
imulations have been done both with non-radiative physics, using 
the same technique and resolutions as the radiative ones. 

We have stored snapshots for 15 different redshifts in the 
range ^ z $C 9 for each resimulated object. In the case of the 
two massive clusters, we stored 320 snapshots for a more detailed 
study of the formation process. A graphical representation of the 
gas distribution for all the MUSIC-2 clusters for each of the 15 
redshifts in radiative and non-radiative runs can be seen online at 
|http : / /music . ft . uam . es| In Fig[T]we show, as an example, 
the gas density, colour coded according to the temperature, and the 
stellar component for 4 different clusters of the MUSIC database at 
z = 0. 

The cumulative mass function of MUSIC-2 clusters in the red- 
shift range ^ z ^ 1 is shown in figure [2] The mass values of the 
dataset completeness for all redshifts analyzed are shown in Table 
[T] spanning from 4.5 x 10 14 /i _1 M at z=\ to 8.5 x 10 14 /i _1 M at 

2 = 0. 

In this first paper about MUSIC clusters, we focused our anal- 
ysis only on the most massive objects of MUSIC-2, selecting all the 
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z 


0.00 


0.11 


0.25 


0.33 


0.43 


0.67 


1.00 


M„ 


8.5 


7.5 


7.5 


7.0 


7.0 


4.5 


4.5 


N 


172 


147 


90 


76 


41 


57 


14 



Table 1. Mass completeness of MUSIC-2 dataset (compared to dark-matter 
only MD simulation) in the redshift range ^ z ^ 1 . The virial mass is re- 
ported in 10 14 h^ 1 Mq and N is the number of clusters beyond the com- 
pleteness mass limit. 




Figure 1. Gas distribution in MUSIC clusters: the physics at z = 0; 
a bullet-like cluster(bottom-left) and a relaxed cluster (bottom-right) of 
MUSIC- 1 at 2 = 0.3. The images of the all MUSIC dataset (available 
at http : //musi c . ft ■ uam.es) have been generated with SPLOTCH 
(Dolag et al.pXjg) . 



clusters of the dataset, both radiative (CSF) and non-radiative (NR), 
with virial mass M v >5xl0 14 /i _1 M©at z — 0. This selection 
allows us to analyze a subset in a mass range where MUSIC consti- 
tutes a complete volume limited sample (more than 80 per cent of 
the clusters found in the 1 h~ Gpc MD box above this mass limit 
have been resimulated by MUSIC). Moreover, those clusters are the 
best resolved ones , containing millions of particles. This allows us 
to to extend the analysis towards the inner regions of the clusters. 
The selected subset has been analyzed at seven different redshifts. 
In Table[T| we show the actual redshifts analyzed and the total num- 
ber of clusters used in each redshift bin. We did not include any of 
the MUSIC- 1 clusters in this analysis, as only very few of them 
have a virial mass larger than 5x 10 14 /i -1 M0at z — 0. 



2.3 Halo-finding and general properties of the clusters 

Haloes and sub-haloes of MUSIC simulations have been iden- 
tified and their properties have been measured using the hybrid 
MPI+OpenMP parallel halo finder AHF (see Knollma nn & Knebe] 
|2009| for a detailed description), a code whose algorithm automati- 
cally identifies haloes, sub-haloes, sub-subhaloes, etc. 

It has been proven that it can reliably recover substructures 
containing at least 20 particles ( Kneb e et al.|201 1\ . 




10 14 
M v [h" 1 M ] 



Figure 2. Cumulative mass function of the MUSIC-2 dataset from z = to 
2=1. The vertical lines indicate the mass limit beyond which the dataset 
constitutes a complete volume limited sample 



The luminosity of the stars formed in the simulation has been 
computed using the STARDUST l |Devriendt et al.||1999| > code, a 
stellar population synthesis model of the spectral energy distribu- 
tions (SEDs) of star bursts from far— UV to radio wavelengths. 
From the SEDs of the stars we could compute the total luminos- 
ity in different photometric bands for each object of the simulation. 

Simulated X-rays maps of our clusters will be produced using 
X-MAS ( |Gardini et al.|2004||Rasia et al.|2008) and PHOX ( fBiffil 
|et al.|2012) , tools developed in order to simulate X-ray observations 
of galaxy clusters obtained from hydrodynamical JV-body simula- 
tions. One of the main characteristics of these codes is to generate 
event files following the same standards used for real observations, 
allowing to analyze simulated observations with the same tools of 
real observations. 

Here, we provide the estimated X-ray temperatures of the sim- 
ulated clusters as: 



Ej T i m i 

J2i mipiA(Ti)Ti 
J2i rmpiAiTi) 

22 i PiTTliTi 



(1) 



(2) 



(3) 



2^i p 1 m l T i 

calculating the mass weighted temperature (Eq[TJ, the emission 
weighted tem perature (Eq|2]> a nd the spectroscopic-like temper- 
ature (Eq[3]( (Mazzotta et al. 2004), where m,:, pi, Ti are the 
gas particle mass, density and electronic temperature and A(T) is 
the cooling function. The particle gas densities of the CSF runs 
have been corrected in order to take into account the multiphase 
model adopted in the radiative simulation; the spectroscopic-like 
and emission-weighted temperatures have been calculated only by 
considering gas particles with kT >0.5 keV and we have assumed 
A(T) oc vT corresponding to thermal bremsstralung. For all the 
MUSIC-2 objects (clusters and groups) the three temperature defi- 
nitions, estimated at the virial radius, are shown in Fig.[3]along with 
the corresponding virial mass, M v , at z = for both the simulations 
(NR and CSF). We will have a more detailed analysis of the cluster 
X-ray properties in a following papers. 
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baryon, gas and star fractions as: 



m. [ir'M s ] 



Figure 3. Distribution of mass-weighted (top panel), emission-weighted 
(middle panel) and spectroscopic-like (bottom panel) temperatures of all 
MUSIC-2 objects (clusters and groups) at virial radius as a function of the 
virial mass at z = 0. Red diamonds refer to CSF objects, blue triangles to 
NR objects. The dashed vertical line indicates the subset analyzed in this 
paper (M v >5x 10 14 /i _1 M at z = 0). 



3 BARYON CONTENT OF MUSIC-2 MASSIVE 
CLUSTERS 

We study the internal baryon content of the MUSIC-2 cluster sam- 
ple estimating the gas and stellar fractions, as well as the total 
baryon fraction. In this way we are able to check the consistency 
between our simulated dataset and available observational results. 
Moreover, as examined in the following section, the gas fraction 
plays a relevant role on the Y — M scaling relation. The over cool- 
ing effect, a hot topic in all cluster simulations (Kravtso v et al.| 
|2005| >, is also monitored in this dataset. 

The mass of each component is estimated inside fixed over- 
density radii defined either with respect to the critical or back- 
ground density. In both cases, we calculate the gas mass M gaa ,A, 
the star mass M s tar,A and the total mass AIa inside a fixed over- 
density radius ta selecting and summing all the particles of gas, 
star and dark matter falling inside that radius. We then define the 



fgas 
J star 
fbar 



M A 

M s tar,A 

Z M A 

M gas> A + M s tar,A 
AfA 



(4) 

(5) 
(6) 



These quantities are calculated for the CSF clusters. In the 
case of the NR simulations, where all the SPH particles are gas, 
only the baryon fraction, defined as fbar=M ga s,A/M& is esti- 
mated. 

From now on, we are going to compare the integrated cluster 
properties at different redshifts. It is therefore very important the 
definition of the aperture radius where the properties are estimated. 
The most often used, more favourable in terms of simplicity, con- 
sists on defining the radius at which the mean enclosed total density 
is a fixed factor (independent on redshift) of the critical density of 
the Universe, p c (z). This is useful to make comparisons with the- 
oretical predictions from the spherical collapse model, which de- 
scribe the virialized part of a cluster in terms of the density contrast 
A„. But A„ varies with redshift and depends on the cosmolog- 
ical parameters adopted. Thus, as we show in Appendix A, this 
definition of aperture radius based on the same critical overden- 
sity value at different redshift does not allow to consider a region 
of the clusters corresponding to the same fraction of the virial ra- 
dius, R v (z) at different redshift. This could introduce a bias in the 
study of the redshift evolution of integrated cluster quantities sim- 
ply due to the different cluster regions probed at different redshifts. 
To alleviate this problem, it is useful to adopt another definition of 
aperture radius based on the redshift-dependent background den- 
sity, pb(z). This approach was introduced the first time by Maughan 
|et al.| {2006} to study the evolution of cluster X-ray scaling rela- 
tions, but it has never been considered in the study of the SZ rela- 
tions before. In this paper we use both definitions of aperture radii: 
a fixed value, independent on redshift, critical overdensity and a 
redshift-dependent background overdensity. We refer the reader to 
Appendix A in which we compare the redshift evolution of the radii 
for these two definitions. 

In order to have a detailed description of the integrated clus- 
ter properties starting from the inner regions up to the virial ra- 
dius, we analyze our clusters at seven different critical overdensi- 
ties, A c =[A„, c (z), 200, 500, 1000, 1500, 2000, 2500], where only 
the virial value is different for each redshift as in Eq.l |A6[ >. Simi- 
larly we select a set of (redshift dependent) background overdensi- 
ties, A b (z)=[A v<b (z), 500(2), 1000(z), 1500(z), 2500(z) , 5000(z), 
7000(2)], where we are adopting the notation used by Maughan 
|et al.|p00"6] >. We are referring to an overdensity At(z) correspond- 
ing to the value at z = but that is changing its value along the 
redshift (i.e. writing 7000(z) it means that we assume A(,=7000 
at z = and varying the value with redshift as in Eq. ^A8[ >). These 
values have been chosen aiming at having the same fraction of the 
clusters volume than those analyzed with a fixed critical overden- 
sity criteria. In this way we can infer the properties at overdensities 
corresponding to cluster regions that are well fitted by X-ray obser- 
vations or by current SZ observations. Due to the high resolution 
of our simulations, we can extend the analysis down to the cluster 
innermost regions (A c = 2500). Overdensities higher than 2500, 
corresponding to the core region of the cluster, can be affected in a 
non negligible way by cold flows or galaxy feedbacks. 

Fig|4] shows the evolution of the gas, star and baryonic 
fractions in a region close to the cluster centre, (A;,=7000(z) 
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Figure 4. Evolution of fgas (red diamonds), fstar (green triangles) and 
fbar (yellow squares) of CSF clusters and fb ar (blue squares) of NR clus- 
ters for A D =7000(z) (top panel) and A c =2500 (bottom panel). 



Figure 5. Evolution of fgas (red diamonds), fstar (green triangles) and 
fbar (yellow squares) of CSF clusters and fb ar of NR clusters (blue 
squares) for A(,=1500(z) (top panel) and A c =500 (bottom panel). 



and A c =2500), corresponding at almost R v /5. When we con- 
sider A c =2500 we observe a moderate decrease in the gas 
fraction as we go from z = (/ 9as =0.09±0.01) to z = 1 
(/ 9as =0.08±0.01). Correspondingly, the value of the star fraction 
grows from / star =0.08±0.01 at z = to / star =0.10±0.01 at z = 1. 
The baryon fraction shows, however, no significant variation with 
redshift, and it is compatible with the cosmic value Qb/Q, m =0.174 
( |Komatsu etal.|2011| > for CSF clusters (/ 6ar =0.174±0.009), while 
it is smaller for NR clusters (/;, ar =0.15±0.01). We find a similar 
scenario when we define the aperture radius based on A(,=7000(,z), 
but in this case the evolution of gas and star fraction with redshift 
is more evident: from / stls =0.09±0.01 and / s t ar =0.07±0.01 at z 
= to / 9as =0.07±0.01 and / star =0.11±0.02 at z = 1. We observe 
that at these overdensities the value of the star fraction is high, es- 
pecially at high redshifts, but the effects of a possible over cooling 
in our simulations do not seem to be strong enough to deviate the 
mean values of the gas fraction from the observational estimates. 

The evolution of the same baryonic component fractions is 
shown in Fig. [5] but for a larger integration domain defined by an 
overdensity: At = 1500(2) and A c = 500 respectively, corre- 
sponding to ~ Rv/2. In this case there are neither any relevant 
differences between the two overdensities used to define the radii, 
nor there is a significant evolution with redshift observed. There- 
fore, we can confidently assume the values at z — for A c =500 
to be representative (within 10% ) of the gas fraction values for 
clusters up to z = 1: / sas =0.118±0.005, / s t a r=0.048±0.003, 
,/w=0.166±0.004. For NR clusters we have / 6ar =0.157±0.006. 

Fig. [6] shows the average radial (normalized to the virial ra- 
dius) cluster profiles of gas, star and baryonic fractions for all the 
redshift considered in this analysis. An evolution of the baryon 
properties with redshift is clearly seen in this plot, which is more 
relevant at smaller R/R v . The gas fraction increases from inner to 
outer regions and shows higher values at low redshifts. On the con- 
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Figure 6. Top panel: Average radial profiles, in units of the virial radius, 
of fbar for CSF clusters (continuous lines) and for NR cluster (dashed) at 
different redshifts from 2=0 (black) to 2=l(red). Bottom panel: f ga s (con- 
tinuous lines) and fstar (dashed lines) profiles for CSF clusters only, from 
z = (black) to z = 1 (red). 
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Figure 7. Comparison among gas fractions of MUSIC-2 clusters at different 
overdensities and observational estimates. Blue dashed area shows the la 
dispersion of f ga s derived from the volume limited NR MUSIC-2 clusters 
while the red solid area corresponds to the CSF MUSIC-2 clusters. The 
x-axis denotes the different observations: Z10 (Zhang et al .|201(ty , V06 
(Vikhlinin et al.|2006> , A08 jAllen et al.|2008) , E10 jEttori et al.|2009) , J10 
(Juett et al.|2010), LR06 , (LaRoque et al.|2006), V09 (Vikhlinin et al.|2009), 
D12jDavidetal.|2012) , M04 jMaughan et al.|2004) , M06 jMaughan et al.| 
|2006| . 



trary, the star fraction decreases when moving towards the virial 
radius and is bigger at higher redshifts, showing that star forma- 
tion is preferentially located at the central parts of clusters at low 
redshift. The baryon fraction decreases when reaching the cluster 
outskirts for CSF clusters (with mean values higher than the cosmic 
ratio) while the NR clusters show the opposite behaviour. In both 
physical configurations the mean baryon fraction becomes smaller 
when moving from z = 1 to z = 0. The evolution of ft ar at large 
fractions of the virial radius is very mild for CSF clusters and con- 
verges towards a nearly constant gas fraction for A c = 500. 

As a summary of the previous discussion we conclude that us- 
ing Ab(jz) to define the aperture radius at different redshift allows 
to reliably compare the same region ( with respect to the virial ra- 
dius of the cluster at each redshift) and thus is the best choice when 
studying redshift evolutions. On the other hand, by selecting the 
radius based on fixed A c values allows us to use a f gas which is 
nearly constant with redshift. 

We now compare the values of the gas fractions derived from 
our simulations with those measured in observations. 

In Fig|7]we show a compilation of published f gas estimates 
from observations at different overdensities and our numerical re- 



sults for the corresponding overdensities, both from the NR and 
CSF MUSIC clusters within the la dispersion values. 

Our results are consistent with the gas fractions that have been 
measured at 3 overdensities (A c = 2500, 1000, 500) by the Lo- 
CuSS project ( |Zhang et al~||2010| >. At A c = 500 we also found 
consistency with f gas estimated by Ettori et al. j2009|,|Juett et al. 



pOlO), [LaRoque et al.|p006},|Vikhlinin et al.| <|2009l>, ~|Ma ughan 
|et al.| ( |2006) , |Maughan et al.| l |2004| > and |David et al.H2012). 

At A c =2500 our results are still consistent with Vikhlinin et al. 
( (2006) and |Maughan et al.| ( [2006) , while | Allen etaT] p008l l shows 
a higher value of the gas fraction. 

At A c =200 our clusters are also compatible with the two 
only measurements at these large radius, Mau ghan et al.| ( |2006| >, 
Maughan et al. (2004), but due to the larger error bars they are less 
restrictive. 



4 SZ CLUSTER SCALING RELATIONS 

The correlation of the SZ signal with other cluster properties was 
initially studied by considering the central comptonization param- 
eter, 2/0, (see |Cooray et al.|1999"l and |McCarthy et al.|2003) but the 
inadequacy of a simple /3-model for the ICM spatial distribution 
was in short clear. For the first time Benson et al. (20041 applied 
the Compton {/-parameter integrated within a redshift-independent 
density contrast, A c =2500, in a SZ scaling relation for 15 clus- 
ters observed with SuZIE II. As the j/o parameter gives only local 
information about the status of the cluster, the integrated Y consti- 
tutes a global property, strongly bound to the energy of the cluster 
l |Mroczkowski et a l. 2009 1 and less model-dependent. 

The comparison between quantities integrated within a con- 
stant critical overdensity is now widely employed in X-ray and SZ 
scaling relations analysis at overdensities A c =2500 (Comis et al 
20TT]|Bonamente et al.|2008] > and A c = 500 jArnaud et al.|2010| 
Planck Collaboration et al. 20 11c Andersson et al. 201 1). The great 
benefit, distinctive in SZ observations of ignoring cluster redshift to 
correct the redshift dependence of the overdensity, is evident. 

In order to give an estimate of the integrated Y of MUSIC 
clusters, we have produced simulated maps of the thermal SZ ef- 
fect, integrating along the projected radius to extract the value of Y 
at fixed redshift and overdensity. 



4.1 Cluster ?/-map generation 

The maps of the Compton {/-parameter are generated with the im- 
plemented procedure employed in |Flores-Cacho et al.H2009) . The 
simulation provides for each i-th SPH particle: position (r{), ve- 
locity density (pi), internal energy (Ui) and SPH smoothing 
length (hi). In CSF simulations we also find the number of ionized 
electrons per hydrogen particle (N e ,i) and metallicity Zc using this 
additional information we could achieve a more detailed estimate 
of electron temperature T e (per gas particle) as: 



T e = (7- l)Um p 



(7) 



where m p is the proton mass, 7 the polytropic coefficient (set to 
5/3 for monatomic gases) and /1 the mean molecular weight of elec- 
trons, that we estimate as 



1 - &He 
1 + 7 + N e 



(8) 
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(Yhe is the nuclear helium concentration). The number electron 
density (in cm -3 ) is derived as 

'1-Z-Yh.* 



(9) 



(The density p gas of the gas particles of CSF particles has been cor- 
rected in order to take into account the multiphase model adopted 
in the simulation) 

We are interested only in the thermal component of the SZ 
effect. This is quantified in the Compton y-parameter: 



V = 



(10) 



The integral can be discretized along the line of sight (los) as 
k B o T 



V ■ 



y^^n Rt jT e! iWp(\ n - r com \,h t ), 



(11) 



where the tSZ signal associated to each particle is spread to the sur- 
rounding area by the mass profile m 9as W p (|ri - r com |, hi), where 
W p is the projection of the normalized spherically symmetric spline 
kernel W used in the simulations: 



W(x,hi 



Tlij 



1 - 6x 2 + 6x 3 , 
2(l~x)\ 



< x sC 0.5 
0.5 < x sC 1 (12) 
x > 1 



4.2 SZE scaling relations 

The integrated Compton y-parameter, Y, appears to be a very good 
cluster mass proxy. With the aim to properly employ SZ observa- 
tions as mass finder, it is fundamental to accurately calibrate the 
Y — M scaling relation. 

As the integrated Compton parameter Y is a very good proxy 
of the clusters mass, it is fundamental to study the Y-M scaling 
relation. The Y parameter is estimated as follows: 



Y = / ydn = D^ 1 



\m e c 2 J J 



dl / n e T e dA 



(13) 



where Da is the angular distance and the integration is performed 
inside the solid angle, fl, equal to an area A in the projected map, 
corresponding to an overdensity A. 

Adopting a self-similar model to describe clusters of galaxies 
(assuming hydrostatic equilibrium and isothermal distribution for 
dark matter and gas particles) the gas temperature is related to the 
mass by a simple power law: 

T e oc (M A E(z)) 2/:i (14) 

where, as we said in the previous section, the total mass is calcu- 
lated inside the overdensity radius rA we are considering and E(z) 
the Hubble constant normalized to its present day value defined as 
E(z)=(n M (l + z)'- i + ^iA + nk(l + z) 2 ) 1/2 . Due to the hypothesis 
of isothermal equilibrium, Y is proportional to the integral of the 
electron density n e over a cylindrical volume, leading to: 



Y C D A <xT e J n e dV oc M gaa T e = f gaa M tot T e (15) 

Considering Eq,^14|>, we find a relation that depends directly on 

Mtot and f gas ; 



YcD\ « f ga sM^E(z) 2/s 



(16) 



To fit a spherical quantity, like the total energy of a cluster, it is 



better to introduce the spherical integrated Ys parameter. This is 
defined as the Compton-y parameter integrated inside a spherical 
volume with a radius equal to rA. The Eq,l|16[> turns out to be: 



Ys,Ac = 



m e c 2 fj, e 



c GHq \ 



2/3 



^(«) 2/3 / 9 a S ,AM 



5/3 



(17) 

where with Ys,a we denote the Y parameter integrated over a clus- 
ter centred sphere of radius r a . Here, we are assuming a constant 
fg as independent on the cluster mass for each overdensity. We cal- 
culate from the map the Ya parameter at the considered overdensity 
integrating along all the pixels of the projected map falling inside 
the projected radius ta- 



Ys,a = ^2 Vi x dl 2 vlx 



(18) 



In order to take into account the different integration domain, we 
define the C parameter as: 



C 



YcDj 
Y s 



(19) 



In Fig[8]we give an estimate of C from the values of Ys and YcD\ 
computed from our clusters at different redshifts and overdensi- 
ties. We find that C is constant with redshift and grows with over- 
density. Available values in literature are consistent with ours: at 
A c =2500 we found C=2.0±0.5, in agreement with the estimate of 
C ~ 2 given b y|Bonamente et al.||2008} and C=1.67±0.80 given 
by |Comis et"aT1<|201 1|>, while a t At;=500 we have C=1.3±0.2, com- 
patible with Andersson et al. ( 201 1 ). From now on we will denote 
the spherical Ys parameter simply as Y, 

The Y — M scaling relation at a fixed overdensity is studied 
performing a best fit of 



Y A = 10" 



M A 
ft- 1 M Q 



Et 



\2/3 



- 2 Mpc 2 } 



(20) 



where Ma is the total mass calculated inside the sphere of radius 
t'a that we are considering and B the normalization, defined as: 



B = 



OT A* 
m e c 2 /x e 



fKcGHo 



2/3 



+ fg. 



(21) 



contains all the constant terms of the first part of the right hand 
of Eq.l|17[l. We remind that, if the clusters follow the self-similar 
model, then the slope A is expected to be 5/3. 

We are also interested in the inverse scaling relation M — Y, 
which is more easily applicable to observational clusters (the mass 
should be inferred from Y measurements): 



Ma = 10 s * 



Ya 



h- 2 Mpc 2 



£(2)- a/B [ft _1 M e 



(22) 



We study these scaling relations in the log-log space, in the 
form logT = B + AlogX, where in the case of the Y — M 
relation we have X — M and T = Y E(z) K f^: k is fixed to its 
expected self-similar scaling with z (k = —2/3) and in this case 
fgas is not included in the normalization parameter B. 

Figure [9] shows the Y — M scaling relation of our dataset at 
A c = 500 for both physical configurations. We choose to consider 
a value of f gas constant at fixed overdensity and redshift, assigning 
to all clusters the mean value of f gaa as: 



fgas(A,z)=J2 



gas,i ^) 



(23) 
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Figure 8. Behavior of the C=Yq D 2 a IY$ parameter with critical overden- 
sity: each value corresponds to the mean value calculated over all the red- 
shifts analyzed. The error bars represent the mean error calculated from the 
mean values of C over all redshifts. 
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Figure 10. The best-fit slope A (top panels) and normalization B (bottom 
panels) of the Y — M scaling relation as a function of overdensity for 2 = 
CSF (left panel) and NR (right panel) MUSIC-2 clusters. 



ied in section[6] we can therefore express the Y—M scaling relation 
at A c = 500 as: 



MUSIC-2(CSF) 
MUSIC-2(NR) 
PLANCK 



' CSFR : A=1.66±0.02 

NR : A=1.67±0.02 
' PLANCK : A=1.74±0.08 
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Figure 9. Y — M scaling relation at z = 0, A c = 500. CSF clusters are 
shown as red diamonds (the best fit as a black solid line) and NR as blue 
triangles (best fit as an orange solid line). The green line shows the best fit 
from Planck data after the removal of the positive bias in the X-ray mass 
estimates with respect to weak lensing masses (see text). 



where N is the number of clusters in the considered subset. 

CSF clusters show higher values than NR clusters, as the rela- 
tion depends on the inverse of f ga s, which is bigger for NR clusters, 
(see Section[3](. 

The values of A and B at z = with their l-a errors are listed 
in Tableland shown in Fig[l0] At all the overdensities considered, 
the values of the slope A seem to agree, within their errors, with the 
self-similar prediction of 5/3. NR clusters seem to show a lower 
slope at high overdensities and CSF clusters an higher one, match- 
ing both with the same self-similar value at A c =500. The value of 
B seems to stay almost constant at all overdensities. 

Leaving aside a possible variation with redshift, which is stud- 



10" 



M 5 



E(z) 2/3 [h- 2 Mpc 2 } 



(24) 

We estimate the scatter on the Y — M relation by computing the 
sum of the residuals as follows: 



O"log 10 Y 



Eti L lQ g Ya.i - (A log Af A ,i + B)Y 
N-2 



(25) 



As shown in Table [2] the scatter on the Y estimate is about 5 per 
cent at all overdensities. 

In the same way we also calculate the M — Y scaling rela- 
tion. Fig.^^shows the relation at z = for A c =500, for which the 
self-similar model predicts a slope A* = 3/5. Table[2]also lists the 
best fit values of A* at z = 0, which seem to be rather independent 
on overdensity and to agree within the errors with the self-similar 
prediction. In Fig[l2]we show the evolution of A* and B* with 
overdensity. We can then express the best-fit M — Y scaling rela- 
tion at A c = 500 as: 



Ah, 



10 1 



Yr M 



h~ 2 Mpc 2 



E(z) 



-2/5 



[/i _1 M G 



(26) 

Moreover, the M — Y relation seems to show a smaller scat- 
ter (3 per cent) than the Y — M. The difference between the two 
relations is shown in Fig |13| in which we show the distribution of 
the residuals at A c = 500. Both distributions are well fitted by a 
Gaussian with zero mean and different widths. The distributions 
of residuals for the M — Y is narrower (i.e. smaller a ) than the 
corresponding one for Y — M. 

Therefore, we can establish that the MUSIC-2 clusters show 
Y — M and M — Y scaling relations that are well consistent with 
the self-similar model both at high and low overdensities. This re- 
sult confirms, if necessary, the validity of the self-similar model in 
using the integrated Y parameter as a proxy of the mass and, at 
the same time, demonstrates the reliability of the MUSIC cluster 
dataset for the study of scaling relations. No significant discrep- 
ancies can be pointed out between the two relations, Y — M and 
M — Y. Therefore we will focus all further analyzes described in 
the next sections only on the first one. 
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Figure 13. Distribution of the residuals of the Y — M (in log Y, left panels) and M — Y (in log M, right panels) scaling relations at z = for A c = 500. 
The red (CSF) and blue (NR) lines represent the best fit Gaussian curves to the residual distribution. 
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Figure 11. M - Y scaling relation at z = for A c = 500: CSF clusters 
are shown as red diamonds (the best fit as a black line) and NR as blue 
triangles (best fit as an orange line). 




Figure 12. Same as Fig|10|but for the M — Y scaling relation. 



4.3 Comparison with observations 

Several scaling-laws have been derived from different SZ obser- 
vations (Bonamente et al. 2008, Andersson & SPT Collaboration 



2010, Arnaud et al. 2010, Comis et al. 201 1 , Planck Collaboration 
|et al.| [201 lb, Andersson et al. 2011). Here, our numerical results 
are compared with the best fit scaling relation Y 5 oo — M500 based 
on 62 nearby clusters detected at high signal-to-noise in the first 
Planck all-sky data set (Planck Collaboration et al. 201 lc). We have 
taken into account the 22±8 per cent mass overestimation of the X- 
ray Planck cluster masses derived from the comparison with weak 
lensing mass estimates l |Planck C ollaboration et al. 2012). This dis- 
crepancy can be attributed to the errors in the mass concentration 
measurements and to an offset in the cluster centres between X-rays 
and lensing observations. Correcting the Planck scaling relation by 
this factor, we find a very good agreement, at A c = 500, with the 
scaling relation of MUSIC-2 CSF clusters, as can be seen in Fig|9] 
The masses assigned to the MUSIC clusters are of course free 
of systematic errors since they are the true masses computed from 
the simulations. We have check how the numerical SZ scaling rela- 
tion would be affected if we assume the hypothesis of hydrostatic 
equilibrium (HSE), to derive the mass of the clusters, as is usually 
done in observations. The HSE assumption constitutes one of the 
most significant uncertainties in the derivation of the observational 
scaling relations. We calculate the HSE cluster mass at a given over- 
density as: 



kTr fd\np d\nT 



GfirriH \dlnr dlnr 



(27) 



We used the mass weighted temperature profile of the gas using 
Eq.lJTJ to calculate the gas particle temperature. Although the spec- 
troscopic like-temperature would have been more fitted to compare 
with X-ray observations, we do not expect that the conclusions 
could be significantly different. The distribution of the Hydrostatic 
Mass Bias (HMB), quantified as AM=(M 500 ,hse-M50o)/M 5 oo, 
can be seen in Fig.[l4] A systematic underestimation of 25 per cent 
of the true mass by using the hydrostatic mass is found, with a dis- 
persion of g = 0.43. The result is consistent with recent numerical 
simulations that also find that the HSE assumption underestimates 
the true mass up to 20 per cent 1 Rasia et al. 2006 , Kay et al. 



Nagai et al.|2007a||Nagai et al.|2007b[|Piffaretti & Valdarnini 



2007 



2008 



Ameglio et al. 2009; Lau et al.|2009 r . A more detailed analysis of 



the HMB, especially at high overdensity regions where a significant 
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Figure 14. Distribution of the hydrostatic mass bias of Afsoo for CSF clus- 
ters at z = 0. The red curve represents the Gaussian curve which best fits 
the distribution. The Gaussian is centered at AM = —0.25 with a width 
of cr = 0.43. 



contribution from non-thermal pressure is present even in relaxed 
clusters, is under way and it will be the topic of a further paper. 



5 THE EFFECT OF THE GAS FRACTION ON THE Y-M 
SCALING RELATIONS 

In the previous section we assumed a constant gas fraction for all 
clusters used to fit the Y — M scaling law at a given redshift and 
overdensity. This assumption is a common observational approach 
in scaling laws study due to the large uncertainty occurring in sin- 
gle cluster measurements. In numerical simulations we can easily 
estimate the gas fraction of each single object and investigate in 
more detail its impact on the Y — M scaling relation. 

In this section we leave the value of f gas free. Therefore, using 
Eq,(|20[> we can derive the modified scaling relation: 



YAfgas 



10 1 



h-^Mr. 



E(z) 2/3 [h 



- 2 Mpc] 



(28) 



In Table [5] we present the best fit values of the Yf ga \ — M 
relation at z = 0, for overdensities in the range A c = 2500 to A„ iC . 
We note that the slope is steeper and the normalization is larger 
(in absolute value) than the Y — M with fixed f ga s, as can be 
seen in Fig. [15] and the distribution of residuals is wider and less 
regular (Fig|16}. Since we saw before that Y — M agrees very well 
with the self-similar predictions, it follows that the departure from 
self similarity can be attributed to the dependence of f gas on the 
cluster total mass. This effect seems to be independent on redshift 
but very sensitive on the overdensity value: core regions present a 
larger deviation from the self-similar values than halo outskirts. 

We now study the relation between / gas — M in the form: 



fgas 



M Ac 
h-^Mr. 



(29) 



It can be easily deduced that introducing the above mass depen- 
dence fg as into the Y — M scaling relation is equivalent to split- 
ting Eq.|20]into two expressions described by Eq.j28} and Eq.l|29[l, 
whose best-fit parameters are connected by the simple relation: 



A = Ai + A 2 
B = B 1 + B 2 



(30) 
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Figure 15. Yf gaB - M scaling relation at z = 0, A c = 2500 for CSF 
clusters (top panel) and NR clusters (bottom panel). The black line shows 
the best-fit relation, the green line the best-fit relation of the Y — M scaling 
relation (with constant f ga s)- 



Since we are dealing with massive clusters, we expect a linear 
relation (in the log space) between the gas fraction and the total 
mass ( |Vikhlinin et al.|2 009 ). In Fig.[l7]we represent the results for 
the fgas — M from our MUSIC clusters for two different overden- 
sities. As can be seen in this figure and in Table [5] the relation is 
steeper at the central parts of the clusters and it is shallower in the 
outskirts while the scatter also increases towards the centre. To test 
for the effects of cooling and star formation in this and other scal- 
ing relations, we always compare them with the results from the 
non radiative simulations of the same clusters. We see from Fig [T7] 
and Table [3} that a similar dependence on mass is also present in 
NR clusters, although somewhat shallower than in the case of CSF 
clusters. 

We can therefore express the best-fit gas fraction dependence 
on mass at A c = 500 as: 

U.04±0.01 

/ fl a S = 10- 16±01 (T^) 02) 



M500 

As we said, this relation becomes much less steep at high overden- 
sities, where the f ga3 — M scaling relation at A c = 2500 can be 
written as: 

0.17±0.01 



fgas 



10 



-3.5±0.2 



M 25 



LTiOIJ 



(33) 



,/i- 1 M 0y 

In Figs.[T8]and[T9]we represent the evolution with redshift of 
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Table 2. Best fit slope, normalization and scatter of Y — M (A, B, <ri og y) and M — Y (A*, B* t cri og m) scaling relations for CSF and NR clusters. 
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Table 3. Best fit slope, normalization and scatter of Yf g 



M (Al,Bl)and/ 9 



M (A2, B2) scaling relations for CSF and NR clusters. 



50 



40 



30 - 



20 



10 



CSF 




NR 




CSF 




NR 



-0.1 0.0 

a log Y 



0.1 



-0.1 



0.0 
a log Y 



0.1 



-0.1 0.0 

a log f„ 



0.1 



-0.1 0.0 
A log f„ 



0.1 



Figure 16. left panel: the distribution of residuals for Yf gas — M at z = and A c =2500. Right panel: same as left panel but for the fgas 
relation. The red (CSF) and blue (NR) lines represent the Gaussian curves which best fit the residual distribution. 
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Figure 17. The f gas — M scaling relation at z = for A c =2500 (left panels) and A c = 500 (right panels). CSF clusters are shown in red diamonds, NR 
clusters in blue triangles. The values of the slope A 2 are also reported in the label. 
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Figure 18. Evolution of the Ai slope parameter of the Yfg~ as ~ M scal- 
ing relation at A c = 2500 (left panel ) and A c =500 (right panel) for CSF 
clusters (red diamonds) and NR clusters (blue triangles). 
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Figure 19. Same as Fig|18|but for the A2 slope. 



the average slopes, A\ and A2, for CFS and NR clusters at two 
overdensities. In Fig.[20]we also show the dependence of the same 
slopes on overdensity for 2 = As can be seen, there is a clear 
linear relation of A\ and A2 on A c . We can therefore parametrize 
the variation of the slopes of the Yfg^s ~~ M and f gaB — M scaling 
relations expressed in equations[28]and[29]as follows: 

Ai = 1.64- 0.65 x (^r) (34) 



A 2 = 0.002 - 0.7 x (35) 

Moreover, from the values of the best fit parameters listed in 
Table [5] the relation of the slopes for the different scaling relations 
expressed by Eq.l|30[l and Eq.|3TJ, is confirmed by our clusters. 



6 REDSHIFT EVOLUTION OF THE Y - M SCALING 
RELATIONS 

In the previous sections we have focused on the study of the Y — M 
scaling relations from the MUSIC clusters at z = 0. In this section 
we will study the possible dependence of the SZ scaling relations 
with redshift. 

We will employ two distinct approaches in our study: 

• Selecting different complete volume limited samples above a 
given mass threshold for different redshift and repeating the study 
done at z — 0. This will give us the variations of the A and B best 
fit parameters of the Y—M relation with redshift ( Sho rt et al.|2010| 
|Fabjan et al.|20TTj|Kay et al.|2012| >. Following previous studies we 
can parametrize the redshift dependence as: 

A(z) = A (l + z) aA (36) 
B(z) = B (l + z) aB (37) 

where Aq and Bo are the values of the slope at z = and ola,b the 
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Table 4. Best fit parameters a a, aj of the evolution of A and B with redshift at different overdensities and physical proccesses. 
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Figure 20. Left panel: the A\ slope as a function of overdensity at 2 = 
for CSF clusters (red diamonds) and NR clusters (blue triangles). Right 
panel: Same as left panel but for the A2 slope. 



best fit parameters that will describe the possible dependence with 
redshift. 

• An alternative approach, more similar to observational stud- 
ies, consisting of generating a sample of clusters that are carefully 
selected from different redshifts, by keeping the same proportion 
of abundance of clusters above a given mass at each redshift con- 
sidered. From these samples, we look for possible dependence on 
redshift of the form: 



Yf-\E(z) 



-2/3 



BM A (l + z 



(38) 



A„(0) 



(39) 



In this case, the normalization parameter B is transformed as 
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A similar approach has been adopted in Ferramacho & Bla nchard] 
< |20 1 1| > with 438 clusters of the X-ray Galaxy Cluster Database 
(BAX) having redshifts ranging between 0.003 and 1 to obtain a 
large sample of potentially observable clusters in SZ. 

Since we want to study the redshift evolution of the scaling 
relations for the two different definitions of aperture radius based 
on the overdensity criteria introduced in section [5] (and discussed 
in more detail in Appendix A), we have to make some preliminary 
considerations. When we use the redshift dependent background 
overdensity to define the integration area, we have to replace A c in 
Eq. ([T7j of the Y - M relation by: 



B = log ■ 



+iog/ gQ j - — -—n m {z) 



A 6 ,„(0)' 



(40) 

But this introduces a natural evolution with redshift due to the scal- 
ing of At(z) that has to be subtracted from the fits to properly 
measure any intrinsic evolution of A and B. 
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Table 5. Number of MUSIC-2 clusters with M„ >5x 10 14 /i _1 M Q at dif- 
ferent redshifts. 



6.1 Evolution with redshift of the Y — M best fit parameters 

In Fig|2T]we represent the mean values of the best fit A and B 
parameters of the complete volume limited mass samples for CSF 
clusters (upper panels) at the same four different overdensities al- 
ready analyzed in section[3] We also estimated the scatter on <ri og y 
for the same cases, calculated as in Eq,l|25|>, 

We can appreciate that the A slope and the normalization B 
show a clear evolution with redshift for low overdensities while 
they flat off when higher overdensity values are considered. The 
situation is the opposite when we consider higher overdensities: at 
these apertures the slope seems to be flatten off up to redshift z ~ 
0.5 and then increases. 

Table|4]shows the value of the best fit parameters a. a and as 
for the same overdensities analyzed in the previous sections. The 
results support the hypothesis of null evolution of Y — AI scaling 
relations at small overdensities (A c =500, Af,=1500(z)) and a pos- 
sible mild evolution at high overdensities (A c =2500, Af,=7000(z)). 

The total mass of the cluster seems to play a fundamental role 
in the evolution of the scaling relations. The same is also certain for 
the f gas — M scaling relation. The Y — M relation seems therefore 
to satisfy self-similarity when considering massive evolved objects 
(we recall that all the clusters analyzed in this paper have M v > 
5x 10 14 /i _1 at z=0). When smaller mass systems are added to the 
catalogues, the deviations from self-similarity become more appar- 
ent, as it is clearly shown in Table. [4] 

The rms scatter of log 10 Y, <xi og y , is also shown in the lower 
row of the upper panel of Fig|2T| Although it is nearly constant at 
low overdensities (about 4 per cent), its value increases (from 6 per 
cent to 8 per cent) at high overdensities and high redshifts. 

If we look at the same evolution in NR clusters (bottom panels 
of Fig|2 1 1>, we find that, in this case, there is no evidence of any 
change in the slope with redshift, and the values of a a and as 
are extremely low at all overdensities (Table[4]l. Also the scatter is 
small (never exceeding 4 per cent) and does not vary with redshift. 

From Fig[2T]and Table [4] we can also conclude that the evo- 
lution of A and B is not affected by the choice of the overdensity 
(fixed critical or redshift dependent background), since both meth- 
ods lead to the same results (evolution at high overdensities and 
smaller masses for CSF clusters and no significant evolution for 
NR clusters). Therefore, there is no need to change from the widely 
used, and more simple, assumption of a fixed critical overdensity to 
define the aperture radius. 
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Figure 21. Evolution with redshift of the slope (top panels), of the normalization (central panels) and of the scatter (bottom panels). The upper part (first three 
rows) of the figure describe CSF clusters, the lower (last three rows) the NR clusters. Redshift varying background density (blue diamonds) at Aj,=1500(z) 
(left) and A 6 =7000(z) (right) and fixed critical overdensity (red diamonds) at A c =500 (left) and A c =2500 (right). 
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Figure 22. Posterior distributions for the parameter at the 4 overdensities considered: A(,= 
spond to 1 and 2 a intervals for the best fit parameters. 



1500(z), A c =500, Ai,=7000(z), A c =2500. The contours corre- 



6.2 Mixed-z sample results 

As we said above, we have also checked the evolution of the Y—M 
scaling relation on mixed-z datasets to mimic the variety of clusters 
that should populate an observational mass-limited galaxy cluster 
survey. These datasets have been built by selecting subsamples of 
CSF objects from all the available redshifts, with the prescription 
that every member appears in the sample only at one redshift and 
that the sample is populated according to the cluster abundances 
of MUSIC-2 as a function of redshift (see Table[5](. This approach 
allows us to evaluate systematics in the reconstruction of the Y—M 
scaling relation. Moreover, if some evolution or mass dependence 
actually exists, it will also help us to estimate the level at which it 
is detectable in a complete mass-limited cluster sample 
We assume Eq.J38l as our reference scaling law: 



=500 



=2500 



A b = 1500(z) A b = 7000(z) 



log(YE(z)- 2/3 f- a \) = B + AlogM + /31og(l + z) (41) 

The fit has been performed with a standard Monte Carlo Markov 
Chain algorithm based on the Metropolis-Hastings sampling 
scheme, where no additional priors have been applied on the 
slope parameter. The convergence has been assessed through the 
Gelman-Rubin test ( Gelman & Rubin 1992 1 on a set of three simul- 
taneous chains run on the same dataset from independent starting 
points in parameter space. The samples were drawn at the four ref- 
erence overdensities (again A c = 2500, 500 and Af,(z)=7000(z), 
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1.672± 0.028 
-29.77± 0.42 
0.17 ±0.10 



1.627 ± 0.022 
-28.85± 0.32 
0.31±0.14 



1.652± 0.028 
-30.65±0.42 
-0.12± 0.11 



1.656 ±0.023 
-30.61± 0.33 
().00± 0.11 



Table 6. Best fit parameters (slope A, normalization B and redshift evo- 
lution /3) for the Y-M scaling relation on a mixed-z sample of 250 CSF 
clusters. 



1500(z)). Fig. [22] shows the posterior distributions of /3 and the 
ft — A and /3 — B joint likelihoods for the reference overdensities. 
The fit results are still fairly consistent with self similarity and no 
additional redshift scaling in the Y — M scaling relation, except for 
the A c = 2500 dataset, which shows marginal evidence of resid- 
ual scaling. Since the mixed-z cluster sample has been drawn from 
the original simulated dataset randomly, the procedure was iterated 
on 5000 samples to check for the possibility of large fluctuations 
in the fit results (e.g. due to outliers in the data distribution). This 
was not the case, since the best-fit parameters at the 4 overdensities 
lied well within the range of fluctuation allowed by the individual- 
sample posterior distributions. In Table [6] all the 3 best fit param- 
eters are listed with their corresponding 1 a errors. The Y — M 
scaling laws for one realization of a mixed-z sample and only for 
low overdensities, are plotted in Figs.|23|and|24| 
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Figure 24. Mixed-z scaling relation at A(, = 1500(z); the different diamond 
colors refer to different redshifts between z = and z = 1 . 



a large number of massive objects. These objects have been sim- 
ulated using two physical processes: non-radiative (NR, gravita- 
tional heating only) and a radiative model (CSF) including cooling, 
UV photoionisation, multiphase ISM and star formation and su- 
pernovae feedbacks. No AGN feedback has been included in the 
radiative simulations. 

We plan to generate mock catalogues from the MUSIC clus- 
ters for the most common observables such as: X-ray and SZE sur- 
face brightness maps, lensing maps and galaxy luminosity func- 
tions. All these deliverables, together with the initial conditions for 
all the MUSIC simulations, will be made publicly available from 
the website http://music.ft.uam.es 

Here we have presented the first results from MUSIC based 
on the analysis of the baryon content and the scaling relations of 
the thermal SZ effect. 

The analysis of the integrated cluster quantities, such as the 
fractions of the different baryon components, or Y parameter, have 
been done using two definitions of the aperture radius. The most 
commonly used approach assumes that the integration domain is 
defined by a fixed value of the mean overdensity inside the cluster 
region compared to the critical overdensity of the Universe. It has 
been shown (see Appendix A) that this definition introduces a bias 
when any integrated quantity is compared at different redshifts. On 
the contrary, by using a redshift-dependent overdensity compared 
to the background density of the Universe we can rescale the aper- 
ture radius such that the integrated region is always the same vol- 
ume fraction of the total cluster volume defined by the virial radius. 

Our numerical results on the baryon content and the SZ scal- 
ing relations have been explored to check possible dependences on 
overdensity, redshift and cluster physics modelling. We have also 
made a comparison between MUSIC clusters and available obser- 
vational results. The scaling relation between SZ brightness, quan- 
tified with the integrated Compton y-parameter, Y, and the cluster 
total mass, M, confirmed the robustness of the self-similarity as- 
sumption for this class of objects present in MUSIC sample. We 
evaluated the effect on the Y — M of having the cluster gas frac- 
tion, f g as, as a free parameter that can vary with the cluster mass. 
Finally, we looked for a possible redshift evolution of the Y — M 
relation using two different approaches. 

Our main conclusions can be summarized as follows: 



7 SUMMARY AND CONCLUSIONS 

Recent new large scale surveys on the far infrared and millimeter 
bands confirmed the validity of the Sunyaev-Zel'dovich Effect as a 
useful cosmological tool to detect and to study galaxy clusters. To- 
day, one of its most relevant applications is to infer the total clus- 
ter mass using the Y — M scaling relation. Recent observations 
(Planck, SPT and ACT) provide the SZ-measurements of a large 
number of clusters from which the Y — M relation can be reliably 
inferred. 

In this paper we have presented the MUSIC dataset, a large 
sample of gasdynamical resimulated massive clusters, that have 
been selected from the MareNostrum Universe simulation and from 
the MultiDark simulation. Our dataset comprises two different sub- 
samples. The MUSIC- 1 which contains 164 objects that have been 
selected according to their dynamical state (i.e. bullet-like clusters 
and relaxed clusters) and the MUSIC-2 which contains a total of 
more than 2000 galaxy groups and clusters. The MUSIC-2 con- 
stitutes a complete volume limited sample for clusters with virial 
masses above 8 x 10 14 h' 1 M Q at z = 0. To our knowledge, there 
are no other works about resimulated clusters that have shown such 



• The mean value of the gas fraction of MUSIC-2 CSF clusters 
at A c =500 is / ga s=0.1 18±0.005. This value is compatible with ob- 



servational results shown in Maughan et al. (20041, LaRoque et al. 
1 2006l|Ma"ughan et al.|(|2006|>,|Vikhlinin et al.|d2009),|Zhang et al. 
| 20T0) , |Ettori et al.| < [2009^ |Juett et al.| j2010| l, and |David et al 



(2012). The values of the baryon fraction at virial radius, for both 



simulation flavours (CSF and NR), are consistent, within the er- 
rors, with the cosmological ratio flj,/fi m (according to WMAP7 
cosmology). At higher overdensities (i.e. A c =2500), we also find 
agreement for the f ga s with observational results from Maughan 
[etaLl ( f2006l ), |Vikhlinin etaT] ( |2009l) and|Zhang et al.| ( |2010^ and~we 
are marginally consistent with jAllen et al.| ( |2008| > and |David et al.| 
( |2012} . This means that the possible effects of numerical overcom- 
ing do not appear to affect dramatically our simulations, despite the 
fact that we have not included strong AGN feedback in our model. 
Concerning the evolution with redshift of the different baryon com- 
ponents, we find significant differences in the evolution depending 
on the definition of aperture radius used. This is more evident at 
large overdensities, where the considered fraction of virial radius 
increases with redshift when a fixed critical overdensity is used to 
the define it (see Fig. |Al[ l. 
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• The Y — M scaling relation derived from our clusters, assum- 
ing a fixed fgas, agrees very well with the predictions of the self- 
similar model and shows a very low dispersion (oi og Y — 0.04). The 
resulting best fit Y 5 oo — M 5 oo relation can then be expressed as: 

1.66±0.02 



Yr >( 



10" 



Mr, 



DO 



E{z)- 2/3 [h- 2 Mpc 2 



(42) 

The M — Y relation, which is more suitable to infer cluster masses 
from Y measurements, is also compatible with the self-similar 
model (which predicts a slope of 3/5) with an even lower dispersion 
(criogM —0.03). The corresponding best-fit M500 — Y500 relation 
can be expressed as: 

0.59±0.01 



M„ 



(43) 

We have compared the above fits with the recent results on the 

Y — M from the Planck Collaboration. The agreement is very good, 
provided that the masses of Planck clusters are overestimated by 
22 per cent due to the biases of the X-ray mass estimations with 
respect to the lensing measurements. We have tested whether this 
bias could be due to the lack of hydrostatic equilibrium hypothesis. 
By comparing the hydrostatic masses with the true one in our MU- 
SIC samples, we find a negative HMB (an opposite trend than in 
Planck clusters). In our case, we find 25 per cent underestimation 
of the true mass by the hydrostatic mass. 

• The dependence of the gas fraction on the cluster total mass 
has been also studied along the cluster aperture radius. There is a 
linear relation between f gas and total mass which is more evident, 
although with larger scatter, when we approach the cluster core. 
This effect is rather insensitive on the adopted physics in the nu- 
merical simulations: both CSF and NR clusters present a similar 
behaviour. 

• Leaving f gas as a free parameter in the Y — M scaling relation 
leads to a deviation from self-similarity. The modified Yf~ aa — M 
relation and the f gaa — M relation are directly related with the 

Y — M if both of them are expressed as power laws. We confirm 
this hypothesis and provide the fits to the f gas — M for two over- 
densities A c = 2500 and A c = 500 

0.17±0.01 



fgas = 10" 



M, 



fgas 



10" 



Mi 



Jr,oo 



(44) 



(45) 



• We have studied the redshift evolution of the slope, and of the 
normalization, in the Y—M scaling law. We did not find any evolu- 
tion for NR clusters. On the contrary, CSF clusters show a marginal 
deviation of the slope from the self-similar prediction at high over- 
density for redshifts larger than 0.5. This result is confirmed using 
two different methods to study the evolution with redshift. 

• No significant differences have been found when comparing 
results at different redshifts using the two alternative definitions 
of aperture areas based on different overdensity criteria: the stan- 
dard, redshift independent, critical overdensity and the redshift- 
dependent background overdensity show the same results, at least 
for the gas fraction analysis, with the latter appearing to be more 
sensitive to the redshift evolution of the sample properties (Y — M 
slope, gas fraction, etc.). We can therefore conclude that even if the 
redshift-dependent background overdensity leads to results which 
better identify a possible evolution, no significant errors are intro- 
duced in the analysis of the SZ effects when using the standard 



critical overdensity criteria. The same conclusions in the case of 
X-ray analyzes have been raised by Maughan et al. 1 2006| > . 

• The use of radiative physics on galaxy cluster simulations in- 
troduces fundamental improvements respect to non-radiative sim- 
ulations: CSF clusters show better agreement with observations in 
the estimate of the gas fraction as well as in the study of the Y — M 
scaling relation. On the other hand, NR clusters overestimate the 
gas fraction and do not seem to be compatible with Y — M scaling 
relations from observational results. 

In summary, we have shown that the MUSIC dataset is well 
suited for the study of massive cluster properties and provides a rea- 
sonable description of observed objects with similar mass. There- 
fore, the MUSIC clusters can be a good cosmological tool. In up- 
coming papers we are going to extend the analysis to other comple- 
mentary scaling relations in X-ray, lensing and optical in order to 
have a full set of observables for a complete volume limited sample 
of ACDM simulated galaxy clusters. 



APPENDIX A: CLUSTER APERTURE RADIUS 
DEFINITIONS 

A correct analysis of the integrated properties of clusters at similar 
or different redshift, is based on the definition of the aperture ra- 
dius. Ideally the virial radius is a suitable choice, confining almost 
the whole cluster. Unfortunately a such large region is challenging 
to be explored due to the lack of sensitivity of present day SZ or 
X-ray observations. Therefore, a smaller aperture radius has to be 
defined. This quantity has to be physically meaningful and must 
hold at different redshifts. 

The widely applied definition of the aperture radius, both in 
observations and simulations, correspond to the radius of a sphere 
that contains a total mass that equals the mass of a sphere with a 
fixed redshift independent overdensity with respect to the critical 
density defined as 



Pc{z) = 



3HjE(z 
SttG 



(Al) 



where Ho is the present value of the Hubble constant, G is the uni- 
versal gravitational constant. 

An alternative definition of the aperture radius consists in 
defining it as the radius of the sphere that contains a total mass 
equal to the mass of a sphere with an overdensity value with re- 
spect to the cosmic background mass density, pb(z) = p c (z)Q m (z) 
where fl m (z) is the cosmological matter density parameter defined 
asfi m (z) = Q M (1 + z) 3 /E 2 (z). 

The cluster centric radius can then be quantified in terms of 
these two overdensities Aj, jC , which means that the cluster average 
density inside rA b is A times the critical or background density. 
The overdensity radius, linked to the critical (or background) den- 
sity, can therefore be expressed as: 



( 3MaM 
\inp c (z)A c (z) 



= ( 3M Ab {z)Q m (z) 
TAb { A7T Pc (z)A b (z) 



1/3 



1/3 



(A2) 



(A3) 



It is important to check whether a redshift independent A c 
could possibly introduce some biases when comparing clusters at 
different redshifts. Our goal is to keep the ratio r/R v constant with 
redshift so that the same fraction of cluster is considered at different 
redshifts. 
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|Maughan et al.| {2006 1 were the first to address this prob- 
lem. They investigated the influence of the definition of overden- 
sity radius only in the case of X-ray scaling relations with a sample 
of 1 1 high-redshift clusters observed with Chandra and/or XMM- 
Newton. Also recently a study of scaling laws specifically, for X- 
ray observables and cluster properties, has faced the same problem 
l |Bohringer et al.|2012) . 

The dependence of the overdensity radius defined in Eq |A2| on 
redshift and on cosmological parameters can be easily derived: 



TA C OC p c (z) 



-1/3 



OC E 



-2/3 



(z) = (n M (i+2) 3 +n A ) 



-1/3 



(A4) 



Similarly, the radius defined with the background overdensity 
depends on redshift as 



r Ah oc 



p c ( 2 )- 1/3 fi m (z)- 1/3 oc(l + Z )- 1 



(A5) 



As can be seen, the dependence on redshift of the two defined 
radius is somewhat different. From z = to z = 1, the ta c de- 
creases by a factor of 1.4 while the rA b decreases a factor of 2. But, 
we are interested to check how the ratio r a c b /Rv depends on red- 
shift. Thus, we take into account the dependence of the overdensity 
at virial radius R v , A v . An analytical fit to the numerical solution 
of the spherical collapse model was given by Bryan & Norman| 
(1998) in a ACDM cosmology (see also |Eke et al.|(l996^ ). 



Av, c (z) ~ 18tt + 82a; - 39a; 2 
as to the background density, is given by: 

A v>b (z) ~ A v , c (z)/(l + x) 



(A6) 



(A7) 



where x — Q, m (z) — 1. 

With this definition of aperture radius the r& c /R v ratio 
changes with redshift as it is shown in right panel of Fig. |A 1 1 The 
effect increases towards the inner regions of the cluster. An over- 
density radius of 2500 over the critical density corresponds to a 
region of 1/ 5th of the virial radius at z = and changes by more 
than 20 per cent towards z — 1. 

Maug han et al.| (2006) suggested an alternative definition of 
aperture radius to show an almost redshift independent variation of 
r/R v . To this end, they propose to define the aperture radius with 
a redshift dependent overdensity value which is a simple rescaling 
of the overdensity defined at z = as: 



A[z) = A(0) 



A v (z) 
A„(0) 



(A8) 



In this way, the ratio of the comoving densities is ensured to be 
constant along the redshift. We refer to this redshift varying over- 
density as A(z), following the notation adopted by Maughan et al 
<20Qg) . 

By looking at the dependence on redshift, cosmology and on 
the density contrast at the virial radius, a background overdensity, 
including the z-dependence (Eq|A8j, seems to be a more suitable 
choice for the A definition, allowing to compare cluster properties 
within radii corresponding to redshift independent cluster fractions. 

b (z)/Rv, 



Al 



, esti- 



In fact, as can be seen in the left panel of Fig. 
mated from the redshift dependent background overdensity, Ab(z), 
looks almost invariant with redshift ( less than 5 per cent variation). 
Thus, applying this definition, we will be able to study the evolu- 
tion of cluster properties almost free from the bias induced by the 
choice of an aperture radius defined by the same overdensity value 
at all redshifts. 




8000 6000 4000 2000 
A b (0) 



2500 2000 1500 1000 500 
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Figure Al. Behaviour of the r/R v ratio at different redshifts (0 ^ z ^ 1) 
for background overdensity varying with redshift (left panel) and for fixed 
critical overdensity (right panel). 
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